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AN INTERACTIVE COMPUTER PROGRAM FOR SIZING SPACECRAFT 
MOMENTUM STORAGE DEVICES 


Floyd J. Wilcox, Jr. 
Langley Research Center 


SUMMARY 

An interactive computer program has been developed which computes the 
sizing requirements for nongimbled reaction wheels, control moment gyros (CMG) , 
and dual momentum control devices (DMCD) used in Earth-orbiting spacecraft. 

The program accepts as inputs the spacecraft's environmental disturbance 
torques, rotational inertias, maneuver rates, and orbital data. From these ^ 
inputs, wheel weights are calculated for a range of radii and rotational speeds. 
The shape of the momentum wheel may be chosen to be either a hoop, solid 
cylinder, or annular cylinder. The program provides graphic output illustrat- 
ing the trade-off potential between the weight, radius, and wheel speed. A 
number of the intermediate calculations such as the X-, Y-, and Z-axis total 
momentum, the momentum absorption requirements for reaction wheels, CMG' s, 
DMCD's, and basic orbit analysis information are also provided as program 
output. 


INTRODUCTION 


Earth-orbiting spacecraft utilize nongimbled reaction wheels, control 
moment gyros (CMG) , and dual momentum control devices (DMCD) for momentum stor- 
age and control, and the development of accurate momentum wheel sizing require- 
ments is essential for effective spacecraft design. As part of the Langley 
Research Center's Computer-Aided Spacecraft Design effort, an interactive com- 
puter program has been developed to size momentum wheels. The program accepts 
as input the spacecraft's environmental disturbance torques, rotational iner- 
tias, maneuver rates, vehicular orientation, and orbital data. Momentum wheel 
weights are calculated for a range of wheel radii and rotational speeds, and 
are provided as graphical output to illustrate the trade-off potential between 
the weight, radius, and wheel speed. Intermediate calculations _ such as the 
X- y- , and Z-axis total momentum, the momentum absorption requirements for 
reaction wheels, CMG's, DMCD's, and basic orbit analysis information are also 
provided as output. 

The momentum wheel sizing computer program listing is given in the 
appendix. 


PROGRAM DESCRIPTION 


The momentum wheel sizing computer program provides an interactive 9 r ^P hi “ 
technique for determining the wheel size and weight for various momentum storag 


& 
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and control devices on Earth-orbiting spacecraft. The program is written in 
ANSI standard FORTRAN on a Control Data Corporation (CDC) CYBER 173 computer. 
The graphics output is generated using Tektronix Plot 10 software routines in 
conjunction with a Tektronix 4014/15 graphic terminal. The program should be 
easily adapted to and can be executed on any host computer with this graphics 
package and remote terminal. 

A simplified flow diagram of the momentum wheel sizing program is shown in 
figure 1. 


Inputs and Assumptions 

The user provides the following inputs: spacecraft environmental distur- 

bance torques, rotational inertias, maneuver rates, orbit altitudes and incli- 
nation, vehicle orientation, and celestial orientation. These inputs are 
provided as outputs from two other Langley Research Center Computer-Aided 
Spacecraft Design programs, (1) the Spacecraft Design and Cost Model (SDCM) , 
reference 1, and (2) the Large Area Space System (LASS), reference 2. The 
orbital equations of motion were obtained from reference 3. The relationship 
between the orbital reference frames and the vehicle reference frames is shown 
in figure 2. The input parameters are: 

INOSE - Vehicle Pointing Orientation 

1. Sideways 

2. Nose Down 

3. Nose Forward 

ISATOR - Celestial Orientation 

1. Earth 

2. Solar 

3. Inertial 

XNNN - Number of Gimbled Configurations for CMG 
A - Apoapsis Altitude 
P - Periapsis Altitude 
ORBINC - Orbit Inclination 

TL - Number of Orbits Between Wheel Unloading 
TACCEL - Maneuver Acceleration Time 
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PDOTR < Y > - X-, Y-, Z-Axis Spacecraft Maneuver Rate 
UPSILN - Angular Pivoting for DMCD 

f M 

TA < Y > - X-, Y-, Z-Axis Atmospheric Disturbance Torque 
TG j^Y^j " x_ » Y- ' z_Axis Gravity Gradient Disturbance Torque 


TS i Y j - X- , Y- , Z-Axis Solar Disturbance Torque 




- X-, Y-, Z-Axis Spacecraft Rotational Inertia 


A major assumption of the program is that the X— axis is the spin axis 
while the Y- and Z-axes are transverse. Also, the spacecraft's orbit may be 
either circular or elliptical. Although all the momentum control wheels are 
normally the same size, the X-, Y-, and Z-axes momentum absorption requirements 
are not identical, therefore, it is necessary to determine the wheel weight for 
the axis with the largest momentum absorption requirement. 


Orbit Analysis 

The first step of the orbit analysis portion of the program is to analyze 
the ellipse representing the orbit, as shown in figure 3. 

From the figure, the semimajor axis, a, is calculated by 
(2R+A+P) 


The radius at periapsis, RP, is calculated by 


RP = R + P 


the focus of the ellipse, c, is calculated by 
c = a - RP 
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and the eccentricity, e, is calculated by 


e 


£ 

a 


where R is the radius of the Earth- 

In the case of a circular orbit, the semimajor axis becomes the circle's 
radius which is also the radius at periapsis. 


a = RP 


This forces both the focus and eccentricity to zero 


c = 0 
e = 0 


The orbital period, T, is then calculated using Kepler's 3/2 Law 


2a 3 / 2 


where u is Earth's gravitational constant. 

The mean orbital motion, N, is calculated using 



and the eccentric anomaly, E, by 


1 9 

E = N + e sin (N) + — e z sin (2N) 

In the case of a circular orbit, the last two terms of the preceding equation 
are forced to zero and the eccentric anomaly becomes 


E = N 
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Using the above calculations, the maximum orbital angular rate of change, 0, 
which is at perigee, is determined from 


0 = 


u\[i - 




e cos 


(E)] 


For a circular orbit, the maximum rate of change is the mean orbital motion 
such that 


0 = N 


Total Momentum 

To calculate the total momentum in each of the three spacecraft's body 
axes, three different momentums are taken into account: (1) maneuver momentum, 

(2) environmental disturbance torque momentum, and (3) orientation tracking 
momentum. 

The maneuver momentums, HXMAN, HYMAN, and HZMAN, in all three axes are 
calculated using the following equations : 


HXMAN 


XJ(PDOTRX) 

57.3 


HYMAN 


YJ(PDOTRY) 

57.3 


HZMAN 


ZJ(PDOTRZ) 

57.3 


To compute the disturbance torque momentum, the total disturbance torques, 
XDT , YDT, and ZDT , are first calculated by summing the atmospheric, solar, 
and gravity gradient torques by 


XDT = TAX + TSX + TGX 


YDT = TAY + TSY + TGY 


ZDT = TAZ + TSZ + TGZ 
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Using these totals, the disturbance torques, HTX, HTY, and HTZ, are then 
computed using 


HTX = XDT (TL) (T) 


HTY = YDT (TL) (T) 


HTZ = ZDT (TL) (T) 


To compute the orientation tracking momentums, HTRAKX , HTRAKY, and 
HTRAKZ , the program first determines whether the spacecraft is solar oriented, 
inertial oriented, or Earth oriented. If the spacecraft is either solar or 
inertially oriented, there is no orientation tracking momentum and the 
following assignments are made: 


HTRAKX = 0.0 


HTRAKY = 0.0 


HTRAKZ = 0.0 


If the spacecraft is Earth oriented, the program will determine whether its 
body axis is sideways, nose down, or nose forward with respect to the orbital 
velocity vector (fig. 2). The following assignments are made if its body axis 

is 


(1) Sideways: 


(2) Nose down: 


(3) Nose forward: 


HTRAKX = XJ (0) 
HTRAKY = 0.0 
HTRAKZ = 0.0 

HTRAKX = 0.0. 
HTRAKY = YJ (0) 
HTRAKZ = 0.0 

HTRAKZ = 0.0. 
HTRAKY = YJ(6) 
HTRAKZ =0.0 


The total momentum in each axis is computed by summing. the maneuver, 
disturbance torque, and orientation tracking momentums 
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HX = HXMAN + HTX + HTRAKX 


HY = HYMAN + HTY + HTRAKY 


HZ = HZMAN + HTZ + HTRAKZ 


Control Device's Momentum Absorption Requirements 


The momentum absorption requirements for the reaction wheel, CMG, and 
DMCD are calculated to determine the wheel weight. The maximum momentum in 
any axis, HMAX, is used for the calculation of nongimbled reaction wheel 
weight. Using an intrinsic FORTRAN function, HMAX is determined by 


HMAX = AMAX1 (HX, HY, HZ) 


Here HMAX will be assigned the largest value between HX, HY, and HZ. To 
compute the CMG momentum for computing wheel weight, the minimum momentum in 
any axis, HMIN, is first calculated using another intrinsic FORTRAN function 


HMIN = AMIN1 (HX, HY, HZ) 


Here HMIN will be assigned the smallest value between HX, HY, and HZ. 
From this, the CMG slew angle, GAMMA, is calculated from 


GAMMA 


TAN 


■i n 


(XNNN-2) HMIN 


(XNNN) HMAX 


Then using the preceding HMIN and GAMMA calculations, the CMG wheel 
momentum, HCMG, is calculated using 


HCMG 


HMIN 

(XNNN) sin (GAMMA) 


The peak gimble rate, DELDOT, is computed by 


DELDOT = 


HCMG 

57.3 (TACCEL) 
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and the peak torquer torque, TCMG, by 


TCMG = 


(HCMG) AMAXl (PDOTRX, PDOTRY , PDOTRZ) 
57.3 


To compute the DMCD momentum for computing wheel weight, the spin and 
transverse axis momentum absorption requirements must be calculated. The DMCD 
configuration for both the spin and transverse axes is shown in figure 4. 

The spin axis absorption momentum, DELTHU, is calculated by 


DELTHU 


HX 

2 


and the transverse axis absorption momentum by 


HU 


AMAXl (HY, HZ) 
^/UPSILN\ 
2 {— 2 — ) 


From these calculations, the total DMCD wheel momentum, HTDMCD , is computed 
by summing the spin and transverse axes absorption momentums 


HTDMCD = DELTHU + HU 


Wheel Weight 


The final step of the program is to calculate the wheel weight. The 
wheel mass is calculated from the relationship between the rotational inertia 
of the wheel and its radius. 


I = MR 2 


Since the angular momentum, L, is equal to the rotational inertia of the 
wheel, I, multiplied by its angular velocity (wheel speed), 


L = IW 
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Then 


L = MR CO 


and 


M = 


r 2 co 


For a solid cylinder: 


I = 


MR 


and 


MR 

L = — — 0) 
2 


and 


M = 


2L 

2 

R CO 


Finally, for an annular cylinder 


I = M 


K 2 ♦ r 2 2 ) 


and 


L = M 


( R l 2 + R 2 2 ) 


CO 


and 


M = 


2L 


CO 


(Rl 2 + R 2 2 ) 



where is the wheel* s inner radius and R 2 is the outer radius. The 

acceleration of gravity/ ACC, at the spacecraft's altitude is then calculated 
by 


ACC = 

RP 


The wheel masses are next multiplied by the acceleration of gravity to deter- 
mine the wheel weight at the spacecraft's altitude. Thus, the equations used 
in computing wheel weight are 


(1) For a hoop 


W = 


L (ACC) 
2 

R 0) 


(2) For a solid cylinder 


W = 


2L (ACC) 
2 

R U) 


(3) For an annular cylinder W 


2L (ACC) 
U)( R 1 2 + R 2 2 ) 


Trade-Off and Output 

Calculated wheel weight curves are shown in figures 5 through 7 for wheel 
radii varying from 1 inch (2.54 cm) to 10 inches (25.4 cm) and wheel speeds of 
500 to 5000 rpm's in increments of 500 rpm's to illustrate the trade-off 
potential between weight, radius, and wheel speed. Weights of hoop, solid 
cylinder, and annular cylinder momentum wheels are presented for the nongimbled, 
CMG and DMCD systems. Calculations were made in U.S. Customary Units. Con- 
version factors for values used in this report are given in table 1. The 
program's input default values, which were used to generate the weight curves 
(figs. 5 to 7) , are given in figure 8. 

In addition to the weight curves, values of 22 calculated parameters are 
outputted, including the maximum torque in any axis which is calculated by 


TMAX 


AMAX1 (HXMAN, HYMAN, HZMAN) 
TACCEL + AMAX1 (XDT, YDT , ZDT) 
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These other calculated parameters are 
RP - Radius at periapsis 
0 - Maximum orbital rate of change 
T - Orbital period 


H 



X-, Y-, Z-axis maneuver momentum 


HT 



X-, Y- , Z-axis disturbance torque momentum 


f x 

HTRAK l Y 

lz 


- x-, Y— , Z-axis orientation tracking momentum 


f X 

H < Y 

lz 


- X-, Y— / Z-axis total momentum 


HMAX - Maximum momentum, any axis 


HCMG - CMG wheel momentum 


HTDMCD - DMCD wheel momentum 


DELDOT - Peak gimble rate 


TCMG - Peak torquer torque 


TMAX - Maximum torque in any axis 


Calculated values of the output parameters are given in figure 9 for the 
weight curves illustrated. 


11 



CONCLUDING REMARKS 


An interactive computer program has been developed at the Langley 
Research Center to size momentum wheels for various momentum storage and 
control devices for orbiting spacecraft. The program considers hoop, solid 
cylinder, and annular cylinder wheels. 

Wheel weights are calculated and are shown for a series of wheel radii 
and wheel rotational speeds. Intermediate calculated parameters are also 
presented. 
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MOMENTUM WHEEL SIZING PROGRAM 

The momentum wheel sizing program listing is given in this 
appendix. 

PROGRAM WHEEL( INPUT, OUTPUT, TAPE1 ,TAPE2,TAPE3) 

REAL N 

COMMON INPUT, LOC, REAL, EXP 

COMMON /WEIR/ WEIGH1 ,WEIGU2,WEIGU3,WEIGH4,WEIGH5, 

*WEIGH6 .WEIGH7 , WEIGHS, WEIGII9 ,WEIGH0 ,R 
COMMON /URPI/U,RP,PI 

DIMENSION R(122), WEIGH! (122) ,WEIGH2(122) 

DIMENSION WEIGH3(122) ,WEIGH4(!22) ,WEIGH5(122) 

DIMENSION WEIGH6(122) ,WEIGII7(122) ,WEIGI 10(122) 

DIMENSION WEIGH9(122) ,WEIGH0(122) 

DIMENSION L0C(5),REAL(9),EXP(12),INPUT(277) 

C RADIUS OF EARTH 

DATA RE/3441.66/ 

C EARTH'S GRAVITATIONAL CONSTANT 

DATA U/1.407850464E16/ 

DATA PI/3.141592654/ 

DATA WEIGHl (1 ) ,UEIGH2(1 ) ,WEIGH3(1 ) ,WEIGH4(1 ) ,WEIGH5(1 ) , 
*WEIGH6(1 ) ,WEIGH7(1 ) ,WEIGH0(1 ) ,WEIGH9(1 ) ,WEIGH0(1 ) , 

*R(1 )/ll*101 ./ 

CALL II1ITT ( 1 20) 

CALL FEET(RE) 

20 PRINT *, "ENTER 1 FOR DEFAULT VALUES" 

PRINT *, "ENTER 2 FOR VALUES FROM LAST RUN" 

PRINT *, "ENTER 3 FOR PERMANENTLY SAVED VALUES" 

PRINT *, "ENTER 4 TO STOP" 

READ *, I VALUE 
PRINT 5025 

IF (IVALUE.EQ.l ) GOTO 40 
IF (IVALUE.EQ.2) GOTO 60 
IF (IVALUE.EQ.3) GOTO 70 
IF (IVALUE.EQ.4) GOTO 999 
GOTO 20 

40 READ (1,5000) (LOC(J) ,J=1 ,3) 

READ (1,5005) (REAL(J) ,0=1 ,9) 

READ (1,5010) (EXP(J) ,J=1 ,12) 

READ (1,5015) (IHPUT(J) ,J=1 ,277) 

REWIND 1 
CALL REPLACE 
GOTO 30 

60 READ (2,5000) (LOC(J) ,J=1 ,3) 

READ (2,5005) (REAL(J) ,J=1 ,9) 

READ (2,5010) (EXP(J) ,J=1 ,12) 

READ (2,5015) (INPUT(J) ,J=1 ,277) 

REWIND 2 
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70 


00 


90 


95 


CALL REPLACE 
GOTO 80 

READ (3,5000) (LOC(U) ,J=1 ,3) 

READ (3,5005) (REAL(J) ,J=1 ,9) 
READ (3,5010) (EXP(J) ,J=1 ,12) 
READ (3,5015) (INPUT(J) ,J=1 ,277) 


REWIND 3 
CALL REPLACE 
INOSE= LOC(l) 
ISATOR= L0C(2) 
XNNN= L0C(3) 
AP0A= REAL(l) 
PER= REAL (2) 
0RBINC= REAL (3) 
TL= REAL(4) 
TACCEL= REAL(5) 
PDOTRX= REAL (6) 
PDOTRY= REAL (7) 
PDOTRZ= REAL(O) 
UPSILN= REAL(9) 
TAX= EXP(l) 

TGX= EXP(2) 

TSX= EXP(3) 
TAY= EXP (4) 

TGY= EXP(5) 
TSY= EXP(6) 
TAZ= EXP (7) 
TGZ= EXP (0) 
TSZ= EXP(9) 


XJ= EXP (10 
YJ= EXP(11 
ZJ= EXP(12) 
WRITE (2,50 
WRITE 


I 

(Ki^OO) (LOC(J).« 
(2,5005) (REAL(J) 


REAL(J) ,0=1 ,9) 

WRITE (2,5010) (EXP(J),J=1,12) 

WRITE (2,5015) (INPUT( J) ,J=1 ,277) 

REWIND 2 

PRINT *,"D0 YOU WANT TO SAVE NEW VALUES PERMANENTLY?" 
PRINT *, "ENTER 1 FOR YES" 

PRINT *, "ENTER 2 FOR NO" 


READ *,JOE 

PRINT 5040 

IF (JOE.EQ.l) 90,95 

WRITE (3,5000) (LOC(J) ,J=1 ,3) 

WRITE (3,5005) (REAL(J) ,J=1 ,9) 

WRITE (3,5010) (EXP(J) ,J=1 ,12) 

WRITE (3,5015) (INPUT(J) ,J=1 ,277) 

REWIND 3 

CALL FEET(APOA) 

CALL FEET(PER) 


14 


o o o 


APPENDIX 


c 

c 

c 

c 

c 

c 

c 

c 


c 


c 


c 

c 

100 


c 

120 

c 

140 


C 

160 


C 


COMPUTE SEMI MAJOR AXIS OF ELIPSE 
A= ( 2 . *RE+APOA+PER) *0 . 5 
COMPUTE ORBITAL PERIOD 
T= 2.*PI*A**(1.5)/SQRT(U) 

COMPUTE RADIUS OF PERIAPSIS 


RP= RE+PER 

COMPUTE FOCUS OF ELIPSE 
C= A-RP 

COMPUTE ECCENTRICITY 


ECC= C/A 

COMPUTE MEAN ORBITAL MOTION 
N= SQRT(U/A**3) 

COMPUTE ECCENTRIC ANOMALLY 
E= N+ECC*SIN(N)+0.5*ECC**2*SIN(2*N) 

COMPUTE MAX ORBITAL ANGULAR RATE 
THETAD= (U*SQRT(1.-ECC**2))/((1.-ECC*C0S(E))**2) 

COMPUTE MOMENTUM REQUIREMENTS FOR SIZING NON-GIMBALLED 


MOMONTUM UHEELS 
MANEUVER MOMENTUM 
HXMAN= XJ*PD0TRX/57.3 
HYMAN- YJ*PD0TRY/57.3 
HZMAN= ZJ*PD0TRZ/57.3 

COMPUTE DISTURBANCE TORQUES 


XDT= TAX+TGX+TSX 
YDT= TAY+TGY+TSY 
ZDT= TAZ+TGZ+TSZ 

COMPUTE DISTURBANCE TORQUE MOMEtITUM 
HTX= XDT*TL*T 
||TY= YDT*TL*T 
UTZ= ZDT*TL*T 

COMPUTE ORIENTATION TRACKING MOMENTUM REQUIREMENTS 
GOTO (120,100,100) ISATOR 

INERTIAL AND SOLAR ORIENTATION 


HTRAKX= 0.0 
HTRAKY= 0.0 
HTRAKZ= 0.0 
GOTO 200 

EARTH ORIENATION 
GOTO (140,160,180) INOSE 
SIDEWAYS 

HTRAKX= XJ*THETAD 
HTRAKY= 0.0 
HTRAKZ= 0.0 
GOTO 200 

NOSE DOWN 
HTRAKX= 0.0 
HTRAKY= YJ*THETAD 
HTRAKZ= 0.0 
GOTO 200 

NOSE FORWARD 
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180 HTRAKX= 0.0 

HTRAKY= YJ*THETAD 
HTRAKZ= 0.0 

C COMPUTE TOTAL MOMENTUM PER ORDIT 

200 HX= (HXMAN+HTX+HTRAKX) 

HY= (HYMAN+HTY+HTRAKY) 

HZ= (HZMAH+HTZ+HTRAKZ) 

C COMPUTE MAX MOMENTUM ANY AXIS 

UMAX= AMAX1 (HX,HY ,HZ) 

C COMPUTE MAX TORQUE ANY AXIS 

TMAX= AMAX1 ( HXMAN , H Yl IAN , HZMAN ) /TACCEL+AMAX1 (XDT.YDT.ZDT) 
C COMPUTE MINIMUM MOMENTUM 

HMIN= AMIN1 (HX,HY,HZ) 

C COMPUTE SKEW ANGLE 

GAMMA= ATAN( (XNNN-2. )*HMIN/(XNNN*HMAX) ) 

C COMPUTE CMG WHEEL MOMENTUfl 

HCMG= HMIN/(XNNN*SIN(GAMMA)) 

C COMPUTE PEAK GIMBLE RATE 

DELDOT= HCMG/TACCEL*57 . 3 
C COMPUTE PEAK TORQUER TORQUE 

TCMG= HCMG*AMAX1 ( PDOTRX , PDOTRY , PDOTRZ ) / 57 . 3 
C COMPUTE SPIN AXIS ABSORPTION REQUIREMENT 

DELTHU= IIX/2. 

C COMPUTE TRANSVERSE AXIS ABSORPTION REQUIREMENT 

HU* AMAX1 (HY,HZ)/(2.*UPSILN/57.3) 

C COMPUTE TOTAL DMCD WHEEL MOMENTUM 

HTDMCD= DELTIIU+HU 

PRINT 5030, (INPUT(J),J=163, 167), APOA 
PRINT 5030, (INPUT(J),J=168, 172), RP 
PRINT 5030, (INPUT(J),J=173, 177), THETAD 
PRINT 5030, (INPUT(J),J=178, 102), T 
PRINT 5025 

PRINT 5035, (INPUT(J),J=273, 277) 

PRINT 5025 

PRINT 5030, (INPUT(J),J=183, 107), HXMAN 
PRINT 5030, (INPUT(J) ,J=188, 192), HYMAN 
PRINT 5030, (IflPUT(J) ,J=1 93, 197), UZMAN 
PRINT 5025 

PRINT 5030, (INPUT(J),J=198, 202), HTX 
PRINT 5030 , ( INPUT(J) ,J=203 ,207 ) ,HTY 
PRINT 5030 , (IflPUT(J) ,J=208,2) 2 ) ,HTZ 
PRINT 5025 

PRINT 5030, (INPUT(J),J=21 3,217), HTRAKX 
PRINT 5030,(INPUT(J),J=218,222),HTRAKY 
PRINT 5030 , ( INPUT( J) , J=223 ,227 ) .HTRAKZ 
PRINT 5025 

PRINT 5030,( INPUT(J) ,J=228,232) ,HX 
PRINT 5030, (INPUT(J),J=233, 237) ,IIY 
PRINT 5030,(INPUT(J) ,J=238,242) ,HZ 
PRINT 5025 
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PRINT 5030, (INPUT(J),J=243, 247), UMAX 
PRINT 5030,( INPUT(J) ,J=248,252) ,TMAX 
PRINT 5030 , ( INPUT ( J ) , J=253 ,257 ) ,HCMG 
PRINT 5025 

PRINT 5030, ( INPUT(J) ,J=258,262) ,DELD0T 
PRINT 5030, (INPUT(J),J=263, 267), TCMG 
PRINT 5030, (INPUT(J),J=268, 272), HTDMCD 
PRINT 5040 
REWIND 3 

220 PRINT *,"D0 YOU WANT WHEEL SPEED, SIZE AND WEIGHT PLOT?" 
PRINT *, "ENTER 1 FOR NON-GIMBLED PLOT" 

PRINT V ENTER 2 FOR CI1G PLOT" 

PRINT *, "ENTER 3 FOR DMCD PLOT" 

PRINT *, "ENTER 4 FOR NO PLOT" 

READ *,IPLOT 

IF (IPLOT.LT.l .0R.IPL0T.GT.4) GOTO 220 
GOTO (240,260,280,20) IPLOT 
240 CALL STUFF ( UMAX, 1 ) 

GOTO 220 

260 CALL STUFF(HCMG,2) 

GOTO 220 

280 CALL STUFF(HTDMCD,3) 

GOTO 220 

999 CALL F 1 1 1 1 TT (0,700) 

5000 FORMAT (515) 

5005 FORMAT (5F10.3) 

5010 FORMAT (5E10.4) 

5015 FORMAT (9(4A10,A5/) ,9(5A10/) ,12(6A10/) ,23(5A10/)) 

5020 FORMAT (5F5.2) 

5025 FORMAT (IX) 

5030 FORMAT (5A10,E14.7) 

5035 FORMAT (6A10) 

5040 F0RI4AT (///) 

STOP 

END 

C 

C SUBROUTINE USED TO CHANGE PLOT AXES VALUES 

C 

SUBROUTINE CHANGE ( XVALUE , YVALUE , XMI N , XMAX , YMI N ,YMAX , J J ) 
INTEGER XVALUE, YVALUE 
5 PRINT *,"D0 YOU WANT NEW GRAPH?" 

PRINT * , " I F YES ENTER 1" 

PRINT * , " I F NO ENTER 2" 

READ *,J0 

IF ( JJ.EQ.l ) 20,10 
10 IF (JJ.EQ.2) 120,5 

20 PRINT *,"00 YOU WANT TO CHANGE AXIS MAX AND- MIN VALUES" 
PRINT *,"IF YES ENTER 1" 

PRINT IF NO ENTER 2" 

READ *,II . 
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IF (II.EQ.l) 40,80 

40 PRINT VENTER X-AXIS MIN AND 11AX VALUE" 

READ *,XMIN,XMAX 

IF ( XMI N . EQ . 0 . 0 . OR . XMI N . GE . XMAX ) GOTO 40 
IF (XMIN.NE.l . 0. OR. XMAX. NE. 10.0) 0J=1 
60 PRINT *, "ENTER Y-AXIS MIN AND MAX VALUE" 

READ *,YMIN,YMAX 
IF (YMIN.GE.YMAX) GOTO 60 

80 PRINT *,"D0 YOU WANT TO CHANGE # OF INTERVALS" 

PRINT *,"IF YES ENTER 1" 

PRINT *,"IF NO ENTER 2" 

READ *,KK 

IF (KK.EQ.l ) 100,120 

100 PRINT *, "ENTER # OF X-AXIS INTERVALS" 

READ *,XVALUE 

PRINT ENTER # OF Y-AXIS INTERVALS" 

READ *,YVALUE 
120 RETURN 
END 
C 

C SUBROUTINE USED TO PLOT CURVES 

C 

SUBROUTINE STUFF (ANGMOM.JJJ) 

INTEGER RPM,XVALUE,YVALUE 

DIMENSION WEIGHl (122) ,HEIGH2(122) ,UEIGH3(122) ,WEIGH4(122) 
DIMENSION WEIGH6(122),WEIGH7(122),UEIGH8(122),UEIGH9(122) 
DIMENSION IJEIGH0(122) ,R(122) ,UE IGM5 ( 1 22 ) ,R1 (122) 

DIMENSION IARRAY(IO) 

COMMON /WEIR/ WEIGH! ,HEIGII2,UEIGH3,UEIGH4,UEIGH5,WEIGH6, 
*WEIGH7, WEIGHS, WEIGH9,UEIGII0,R 
COMMON /URPI/ U,RP,PI 
XMI N« 1.0 
XMAX= 10.0 
YMIN= 0.0 
YMAX= 5.0 
XVALUE= 9 
YVALUE= 5 

20 PRINT *,"D0 YOU WANT WHEEL PLOT FOR:" 

PRINT *," 1- HOOP" 

PRINT 2- SOLID CYLINDER" 

PRINT 3- ANNULAR CYLINDER" 

READ *,M 

IF (M.LT.l .0R.M.GT.3) 20,40 
40 IF (M.EQ.3) 60,80 

60 PRINT *, "ENTER THICKNESS OF ANNULAR CYLINDER IN INCHES" 

READ *, THICK 
80 1=1 

ZNUM= (XMAX-XMIN)/100. 

RADIUS= XHIN 
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IF (M.EQ.3) 100,140 
100 DO 120 J=2,102 
R(J)= RADIUS 
R1(J)= RADIUS+THICK 
RADIUS 3 RADIUS+ZNUM 
120 CONTINUE 
GOTO 180 

140 DO 160 J=2,102 
R(J)= RADIUS 
R1(J)= 0.0 
RADIUS 3 RADIUS+ZNUM 
160 CONTINUE 
180 ACC 3 U/RP**2 

IF (M.EQ.l ) 200,220 
200 HUMACC 3 ANGM0M*ACC 
GOTO 240 

220 HUMACC 3 2.*ANGM0M*ACC 
240 DO 480 RPM= 500,5000,500 
OMEGA 3 RPM/60.*2.*PI 
DO 460 J=2,102 

GOTO (260,280,300,320,340,360,380,400,420,440) I 
260 UEIGH1 (J)= HUMACC/ ( OMEGA* (R(J)**2+R1 (J)**2)) 

GOTO 460 

280 WEIGH2(J)= HUMACC/ ( OMEGA* ( R ( J ) **2+Rl ( J ) **2 ) ) 

GOTO 460 

300 WEIGH3(U)= HUMACC/ ( OMEGA* (R( J)**2+R1 (J)**2) ) 
GOTO 460 

320 WEIGH4(J) 3 HUMACC/ ( OMEGA* ( R( J ) **2+Rl (J)**2)) 
GOTO 460 

340 WEIGH5(J)= HUMACC/ ( OMEGA* ( R ( J ) **2+Rl ( J ) **2 ) ) 
GOTO 460 

360 WEIGH6(J)= HUMACC/ ( OMEGA* ( R( J)**2+R1 (J)**2) ) 
GOTO 460 

380 WEIGH7(J)= HUMACC/ (OMEGA* ( R( J ) **2+Rl (J)**2)) 

400 WEIGH8(J) 3 HUMACC/ ( OMEGA* (R(J)**2+R1 (J)**2)) 
GOTO 460 

420 WEIGH9(J)= HUMACC/ (OMEGA* (R(J)**2+R1 ( J)**2) ) 
GOTO 460 

440 WEIGHO(J) 3 HUMACC/ ( OMEGA* (R(J)**2+R1 (J)**2) ) 

460 CONTINUE 
1= 1+1 

480 CONTINUE 
490 CALL BINITT 

CALL SLIMYO 55,730) 

CALL SLIMX(1 70,920) 

CALL XT ICS (XVALUE) 

CALL YTICS (YVALUE) 

CALL XNEAT (0) 

CALL YNEAT (0) 
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CALL NEWPAG 

CALL NOTX (5,500,5, 5HMHEEL) 

CALL NOTX ( 0 ,485 ,6 ,6HWEIGHT) 

CALL N0TX(8,470,4,4H(LB)) 

CALL NOTX (440,60,17,17HWHEEL RADIUS (IN)) 

CALL NOTX (370,0, 36, 36HRPM'S VARY FROM 500 TO 5000 STEP 500) 
GOTO (500,520,540) JJJ 

500 CALL NOTX (445,30,16,16IIN0N-GIMBLED PLOT) 

GOTO 560 

520 CALL NOTX (480,30,8,8HCMG PLOT) 

GOTO 560 

540 CALL NOTX (470,30,9,9UDMCD PLOT) 

560 GOTO (580,600,620) M 

580 CALL NOTX (490,1 5,4,4HH00P) 

GOTO 640 

600 CALL NOTX (450,15,14,14IIS0LID CYLINDER) 

GOTO 640 

620 CALL NOTX (350,1 5,28,28UANNULAR CYLINDER: THICKNESS-) 

CALL NOTX (650,1 5, 6, CHINCHES) 

CALL FFORM (THICK,5,2,IARRAY,32) 

CALL NOTAT (600,1 5, 5, I ARRAY) 

640 CALL DLIMX (XMIN.XMAX) 

CALL DLIMY (YMIN.YMAX) 

CALL CHECK (R.UEIGH1) 

CALL DSPLAY (R.WEIGH1) 

CALL CPLOT (R.UEIGH2) 

CALL CPLOT (R.WEIGH3) 

CALL CPLOT (R.UEIGH4) 

CALL CPLOT (R, WEIGHS) 

CALL CPLOT (R.UEIGH6) 

CALL CPLOT (R.WEIGH7) 

CALL CPLOT (R, WEIGHS) 

CALL CPLOT (R.WEIGII9) 

CALL CPLOT (R.WEIGHO) 

CALL HOME 
CALL EPAUSE 

CALL CHANGE (XVALUE,YVALUE,XI1IN,XMAX,YMIN,YMAX,JJ) 

IF (JJ.EQ.2) RETURN 
IF (JJ.EQ.3) GOTO 490 
GOTO 80 
END 
C 

C SUBROUTINE USED TO CHANGE NAUTICAL MILES TO FEET 

C 

SUBROUTINE FEET (CHANGE) 

CHANGE® CHAtlGE*6076.103 
RETURN 
END 
C 
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c 

c 


20 


40 


60 


80 


100 


120 


140 

160 


180 


200 


220 

240 


SUBROUTINE USED TO DISPLAY INPUTS 

SUBROUTINE REPLACE 

COMMON INPUT, LOC , REAL , EXP 

DIMENSION L0C(5),REAL(9),EXP(12),INPUT(277) 

PRINT 5000,INPUT(1 ) ,INPUT(2) ,INPUT(3) , 
*INPUT(4),INPUT(5),L0C(1) 

PRINT 5015JNPUT(J) ,INPUT(J+1 ) ,INPUT(J+2) , 
*INPUT(J+3),INPUT(J+4) 

PRI HT 1 ’ 5005 , 1 IIPOT ( 2 1 ) ,IHPUT(22) ,INPUT(23) , 
*INPUT(24) ,INPUT(25) ,L0C(2) 

PRINT 501 5*1 NPUT ( J ) , INPUT ( J+l ) , I NPUT ( J+2 ) , 
*INPUT(J+3),INPUT(J+4) 

PRINT^Ol 0 ,INPUT(41 ) ,INPUT(42) , INPUT (43) , 
*INPUT(44),INPUT(45),L0C(3) 

PRINT 5020 
READ *, CHECK 
PRINT 5035 

IF (CHECK. EQ. 1 ) GOTO 120 

IF (CIIECK.EQ.2) 100,80 

CALL ONE 

PRINT 5035 

GOTO 20 

11= 46 

PRI NT° 5025 *, J , I NPUT ( 1 1 ) , I NPUT ( 1 1 +1 ) , I NPUT ( 1 1+2 ) , 

*INPUT(II+3) ,INPUT(II+4) ,REAL(J) 

11= I 1+5 
CONTINUE 
PRINT 5020 
READ *, CHECK 
PRINT 5035 

IF (CHECK. EQ. 1 ) GOTO 200 

IF (CHECK. EQ. 2) 180,160 

CALL TWO 

PRINT 5035 

GOTO 120 

11= 91 

PRlMT°503o|o^IUPUT(II) »IUPUT(II+1 ) »INPUT(II+2), 

*INPUT( I 1+3) , INPUT ( I 1+4) ,INPUT(II+5) ,EXP(J) 

11= I 1+6 
CONTINUE 
PRINT 5020 
READ *, CHECK 
PRINT 5035 



APPENDIX 


IF ( CHECK. EQ.l) GOTO 280 
IF (CHECK.EQ.2) 260,240 
260 CALL THREE 
PRINT 5035 
GOTO 200 
280 RETURN 

5000 FORMAT (IX, "1 ",1X,4A10,A5,I3) 

5005 FORMAT (1X,"2",1X,4A10,A5,I3) 

5010 FORMAT (1X,"3",1X,4A10,A5,I3) 

5015 FORMAT (3X,4A10,A5) 

5020 FORMAT (///,1X,"IF INPUTS OK ENTER 1",/,1X,"IF WANT TO CHANGE", 
*" ENTER 2") 

5025 FORMAT (I2,1X,5A10,F10.3) 

5030 FORMAT (I2,1X,6A10,E10.4) 

5035 FORMAT (//) 

END 

C 

C SUBROUTINE USED TO CHANGE INTEGER INPUTS 

C 

SUBROUTINE ONE 

COMMON INPUT, LOC, REAL, EXP 

DIMENSION LOC ( 5 ) ,REAL(9) ,EXP(12) ,INPUT(277) 

20 PRINT *, "ENTER NUMBER, NEW VALUE (INTEGER)" 

PRINT *, "ENTER 0,0 TO STOP" 

40 READ *, NUMjIVALUE 
IF (NUM.EQ.O) RETURN 

IF ((NUM.EQ.l .0R.NUM.EQ.2).AND. (IVALUE.LT.l .0R.IVALUE.GT.3)) 
*GOTO 20 

IF (NUM.LT.l .0R.NUM.GT.3) GOTO 20 
LOC ( HUM) = I VALUE 
GOTO 40 
RETURN 
END 
C 

C SUBROUTINE USED TO CHANGE REAL INPUTS 

C 

SUBROUTINE TWO 

COMMON INPUT, LOC, REAL, EXP 

DIMENSION L0C(5) ,REAL(9) ,EXP(12) ,INPUT(277) 

20 PRINT *, "ENTER NUMBER, NEW VALUE (REAL)" 

PRINT *, "ENTER 0,0. TO STOP" 

40 READ *,NUM, VALUE 

IF (NUM.EQ.O) RETURN 
IF ( NUM.LT.l. OR. NUM.GT. 9) GOTO 20 
REAL(NUM)= VALUE 
GOTO 40 
RETURN 
END 
C 
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c 

c 


20 

40 


/ 


SUBROUTINE USED TO CHANGE EXP. INPUTS 


SUBROUTINE THREE 
COMMON INPUT ,LOC,REAL,EXP 
DIMENSION L0C(5),REAL(9),EXP(12),IHPUT 277) 
PRINT ENTER NUMBER, NEW VALUE (EXP)" 
PRINT *, "ENTER 0,0. TO STOP" 

READ *,NUM, VALUE 

IF (NUM.EQ.O) RETURN 

IF (NUM.LT.l .0R.NUM.GT.12) GOTO 20 


EXP(NUM)= VALUE 
GOTO 40 
RETURN 
END 
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TABLE I. CONVERSION FACTORS FROM U.S. CUSTOMARY TO S.I. UNITS 


1 inch 


2.54 centimeters 

1 foot 

= 

0.305 meters 

1 mile 

= 

-1.61 kilometers 

1 pound (wt) 

= 

4.95 newtons = 0.454 kilograms 

1 foot-pound 

= 

1.38 newton-meters 

1 slug-foot 2 

= 

2 

1.36 kilogram-meters 

2 

1 slug-foot 

= 

2 

1.36 kilogram-meters 

seconds 


seconds 





Calculate Weight of Wheel 

— [ 

I 

Plot Wheel Radius versus Wheel 
Weight for Various Wheel Speeds 


Figure 1. Flow diagram. 








NOSE FORWARD 


NOSE DOWN 


SIDEWAYS 


Figure 2 


Relationship between the orbital reference frame and the 
three vehicle reference frames. 







LEGEND: 

a = ^emimajor axis of orbit ellipse 

A = Apoapsis altitude 

C = Focus of orbit ellipse 

P = Periapsis altitude 

R = Radius of Earth 

RP = Radius of periapsis 



Figure 3. Orbit definition. 
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ABSORPTION ABSORPTION 


Figure 4. DMCD configuration. 
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WHEEL RADIUS (IN) 

RPH'S VARY FROM SO TO 5SS» STEP 509 

! (a) Hoop._ 

Figure 5. Weight curves for nongimbled reaction wheel 
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(c) Annular cylinder. 


Figure 5. Concluded 
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(b) Solid cylinder. 


Figure 6. Continued. 
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Figure 6. Concluded 
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Figure 7. Weight curves for DMCD 
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Fisure 7. Continued. 
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Figure 7. Concluded • 


1 INOSE- VEHICLE POINTING ORIENTATION 3 

1= SIDEWAYS 
2= NOSE DOUH 
3= HOSE FORWARD 

2 ISATOR- CELESTIAL ORIENTATION 1 

1= EARTH 
2= SUN 
3= INERTIAL 

3 XNNN- NUMBER OF GIMBALED CONFIGURATIONS 4 


IF INPUTS OK ENTER 1 
IF WANT TO CHANGE ENTER 2 


* 

1 APOA- 

2 PER- 

3 ORBINC- 

4 TL- 

5 TACCEL- 

6 PDOTRX- 

7 PDOTRY- 

8 PDOTRZ- 

9 UPSILN- 


LENGTH OF APOAPSIS (ML) 

LENGTH OF PERIAPSIS (I1L) 

ORBIT INCLINATION (DEG) 

TIME BETWEEN WHEEL UNLOADING (ORBITS) 
MANEUVER ACCELERATION TIME (SEC) 
X-AXIS MANEUVER RATE (DEG/SEC) 

Y-AXIS MANEUVER RATE (DEG/SEC) 

Z-AXIS MANEUVER RATE (DEG/SEC) 
ANGULAR PIVOTING FOR DMCD (DEG) 


175.000 

100.000 

57.300 

1.000 

20.000 

1.000 

1.000 

1.000 

20.000 


IF INPUTS OK ENTER 1 
IF WANT TO CHANGE ENTER 2 


1 TAX- 

X-AXIS 

2 TGX- 

X-AXIS 

3 TSX- 

X-AXIS 

4 TAY- 

Y-AXIS 

5 TGY- 

Y-AXIS 

6 TSY- 

Y-AXIS 

7 TAZ- 

Z-AXIS 

8 TGZ- 

Z-AXIS 

9 TSZ- 

Z-AXIS 

10 XJ- 

X-AXIS 

11 YU- 

Y-AXIS 

12 ZJ- 

Z-AXIS 


ATMOSPHERIC DIST. TORQUE (FT-LBS) 0 

GRAVITY GRAD. DIST. TORQUE (FT-LBS) 

SOLAR DIST. TORQUE (FT-LBS) 0 

ATMOSPHERIC DIST. TORQUE (FT-LBS) 

GRAVITY GRAD. DIST. TORQUE (FT-LBS) 

SOLAR DIST. TORQUE (FT-LBS) 

ATMOSPHERIC DIST. TORQUE (FT-LBS) 

GRAVITY GRAD. DIST. TORQUE (FT-LBS) 

SOLAR DIST. TORQUE (FT-LBS) 

SPACECRAFT ROTATIONAL INERTIA (SLUG-FT**2) 
SPACECRAFT ROTATIONAL INERTIA (SLUG-FT**2) 
SPACECRAFT ROTATIONAL INERTIA (SLUG-FT**2) 


IF INPUTS OK ENTER 1 
IF WANT TO CHANGE ENTER 2 


Figure 8. Input default values. 


!2859E-10 

Is541 E-06 
.3791 E-10 

.8477E-07 


.1595E+04 

.1380E+04 

.8034E+03 
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APOA- APOAPSIS ALTITUDE (FT) 

RP- RADIUS AT PERIAPSIS (FT) 

THETAD- MAXIMUM ORBITAL ANGULAR RATE (RAD/SEC) 
T- ORBITAL PERIOD (SEC) 

ALL MOMENTUM IN (SLUG-FT**2)/SEC 

HXMAN- X-AXIS MANEUVER MOMENTUM 
HYMAN- Y-AXIS MANEUVER MOMENTUM 
HZMAN- Z-AXIS MANEUVER MOMENTUM 

HTX- X-AXIS DISTURBANCE TORQUE MOMENTUM 
HTY- Y-AXIS DISTURBANCE TORQUE MOMENTUM 
HTZ- Z-AXIS DISTURBANCE TORQUE MOMENTUM 

HTRAKX- X-AXIS ORIENTATION TRACKING MOMENTUM 
HTRAKY- Y-AXIS ORIENTATION TRACKING MOMENTUM 
HTRAKZ- Z-AXIS ORIENTATION TRACKING MOMENTUM 

HX- TOTAL X-AXIS MOMENTUM 

NY- TOTAL Y-AXIS MOMENTUM 

HZ- TOTAL Z-AXIS MOMENTUM 

UMAX- MAXIMUM MOMENTUM ANY AXIS 

TMAX- MAXIMUM TORQUE ANY AXIS (FT-LBS) 

hcmg- cmg wheel momentum 

DELDOT- PEAK GIMBAL RATE (RAD/SEC) 

TCMG- PEAK TORQUER TORQUE (FT-LBS) 

HTDMCD- DMCD WHEEL MOMENTUM REQUIR. 


.1 06331 8E+07 
.2151 949E+08 
.1194795E-02 
.5370454E+04 


.2783595E+02 

.2408377E+02 

.1402094E+02 

.1 53541 3E-06 
.45871 08E-02 
.4552533E-03 


0 . 

.1 64881 8E+01 

0 . 


.2783595E+02 
.257371 7E+02 
.1402140E+02 

.2783595E+02 
. 1 391 798E+01 
.1 435261 E+02 

.4112024E+02 
. 250481 9E+00 
.5078648E+02 


Figure 9. Calculated output values. 




1. Report No. 

NASA TM 81798 


4. Title and Subtitle 


2. Government Accession No. 


AN INTERACTIVE COMPUTER PROGRAM FOR SIZING 
SPACECRAFT MOMENTUM STORAGE DEVICES 


7. Author(s) 

Floyd J. Wilcox, Jr.* 


9. Performing Organization Name and Address 

NASA Langley Research Center 
Hampton, VA 23665 


12. Sponsoring Agency Name and Address 

% National Aeronautics and Space Administration 
Washington, DC 20546 

15. Supplementary Notes 


3. Recipient's Catalog No. 


5. Report Date 

March 1980 


6. Performing Organization Code 


8. Performing Organization Report No. 


10. Work Unit No. 
540-02-13 


11. Contract or Grant No. 


13. Type of Report and Period Covered 

Technical Memorandum 

14. Sponsoring Agency Code 


^Cooperative education student at Georgia Institute of Technology 


16. Abstract 


An interactive computer program has been developed which computes the sizing 
requirements for nongimbled reaction wheels, control moment gyros (CMG) , and 
dual momentum control devices (DMCD) used in Earth-orbiting spacecraft. The 
program accepts as inputs the spacecraft's environmental disturbance torques, 
rotational inertias, maneuver rates, and orbital data. From these inputs, wheel 
weights are calculated for a range of radii and rotational speeds. The shape of 
the momentum wheel may be chosen to be either a hoop, solid cylinder, or annular 
cylinder. The program provides graphic output illustrating the trade-off 
potential between the weight, radius, and wheel speed. A number of the inter- 
mediate calculations such as the X-, Y-, and Z-axis total momentum, the momentum 
absorption requirements for reaction wheels, CMG's, DMCD's, and basic orbit 
analysis information are also provided as program output. 


17. Key Words (Suggested by Author(s)) 

Momentum storage 
Spacecraft design 


19. Security Classif. (of this report) 
Unclassified 


18. Distribution Statement 


Unclassified - Unlimited 


20. Security Classif. (of this page) 
Unclassified 


Subject Category 61,18 


21. No. of Pages 22. Price* 

39 $ 4.50 


N-305 


For sale by the National Technical Information Service, Springfield, Virginia 22161 














